Load R packages.

library(Ecdat)
library(HistData)
library(ggvis)
library(plotly)
library(sf)
library(tidyverse)
library(viridis)

3 ways to make maps with R

With ggplot:

ggplot(data=Snow.streets) + 
   geom_path(aes(x=x, y=y, group=street))

With ggvis:

Snow.streets %>%
  group_by(street) %>% 
  ggvis(x=~x, y=~y) %>%
  layer_paths()

With plotly:

plot_ly(Snow.streets, type="scatter", showlegend=F) %>% 
  group_by(street) %>%
  add_lines(x=~x, y=~y)

Improving the plot:

plot_ly(data=Snow.streets, type="scatter") %>%
  group_by(street) %>%
  add_lines(x=~x,
            y=~y,
            hoverinfo="none",
            line=list(color="rgba(0, 0, 0, 1)"),
            showlegend = F) %>%
  add_markers(data=Snow.deaths,
              x=~x,
              y=~y,
              hoverinfo="none",
              marker=list(symbol=0, size=4, color="rgba(153, 0, 0, 1)"),
              name="cholera\ndeath") %>%
  add_markers(data=Snow.pumps,
              x=~x,
              y=~y,
              hoverinfo="text",
              text=~label,
              marker=list(symbol=2, size=10, color="rgba(0, 0, 153, 1)"),
              name="water\npump") %>%
  layout(xaxis=list(visible=F),
         yaxis=list(visible=F),
         title="1854 London")

Making choropleth maps

us_map <- map_data("state")
us_prod <- Produc[Produc$year == 1985, ]
us_prod$region <- gsub("_", " ", tolower(us_prod$state))
merged_data <- left_join(us_map,
                         us_prod[, c("region", "gsp")],
                         by="region")
choropleth <- ggplot(data=merged_data) +
  geom_polygon(aes(x=long, y=lat, group=group, fill=gsp), color="black") +
  theme_void() +
  theme(legend.position="bottom") +
  scale_fill_continuous(guide=guide_colorbar(barheight=unit(2, units="mm"),
                                             barwidth=unit(5, units="cm")))
choropleth

Zooming in on the choropleth map:

tristate <- c("new york", "new jersey", "connecticut")
long_lim <- merged_data$long[merged_data$region %in% tristate]
lat_lim <- merged_data$lat[merged_data$region %in% tristate]
outer <- merged_data[!merged_data$region %in% tristate, ]
choropleth +                                                                            geom_polygon(data=outer,
               aes(x=long, y=lat, group=group),
               fill="gray",
               color = 'black') +
  coord_fixed(xlim=c(min(long_lim), max(long_lim)),                                                 ylim=c(min(lat_lim), max(lat_lim)),                                                   ratio=1.3)

Customizing the map:

h1 <- ggplot(data=merged_data) +                                                  
  geom_polygon(aes(x=long, y=lat, group=group, fill=gsp), color="black") +
  theme_void()
clrs <- rev(viridis::magma(8))
brks <- c(0, 10000, seq(100000, 500000, 100000))
lbls <- format(brks, nsmall=2, big.mark=",", scientific=F)
bar <- guide_colorbar(barheight=unit(4, units="cm"), 
                      barwidth=unit(2, units="mm"))
h2 <- h1 +
  scale_fill_gradientn(colors=clrs,
                       breaks=brks,
                       labels=lbls, 
                       guide=bar,
                       name="GSP (US$)") 
h3 <- h2 +
  labs(title="Gross State Product (GSP) - 1985", 
       subtitle="United States of America") +
  theme(legend.position=c(.9, .1), 
        text=element_text(face="bold"), 
        plot.title=element_text(hjust=.5), 
        plot.subtitle=element_text(hjust=.5))
## Warning: A numeric `legend.position` argument in `theme()` was deprecated in ggplot2
## 3.5.0.
## ℹ Please use the `legend.position.inside` argument of `theme()` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
h3

Reading shape files

districts <- read_sf("../../data/district-boundaries/admin3_poly_32.shp")
health <- read_sf("../../data/health-facilities/all_bphs.shp")    
bg <- element_rect(fill="lightgray")
ggplot() +
  geom_sf(data=districts,
          color="white", fill="lightblue", size=.1) +
  geom_sf(data=health,
          color="red", shape=3) +
  theme_bw() +
  theme(panel.background=bg)

  labs(title="Afganisthan's Health Facilities",
       subtitle="source: www.mapcruzin.com")
## $title
## [1] "Afganisthan's Health Facilities"
## 
## $subtitle
## [1] "source: www.mapcruzin.com"
## 
## attr(,"class")
## [1] "labels"

Interactive globe

data(bankingCrises)
row_year <- match(2009, bankingCrises$year)
countries <- colnames(bankingCrises[, bankingCrises[row_year, ]==1])
axis_props <- list(showgrid=T,
                   gridcolor=toRGB("gray40"),
                   gridwidth = 0.5)
globe <- list(
  showland=T,
  showlakes=T,
  showcountries=T,
  showocean=T,
  countrywidth=0.5,
  landcolor=toRGB("grey90"),
  lakecolor=toRGB("white"),
  oceancolor=toRGB("white"),
  projection=list(type="orthographic",
                  rotation=list(lon=-100, lat=40, roll=0),
                  lonaxis=axis_props,
                  lataxis=axis_props)
  )
plot_geo(width=528, height=528, locationmode="country names") %>%
  add_trace(locations=~countries,
            showscale=F,
            z=1,
            hoverinfo="text",
            text=~countries) %>%
  layout(geo=globe)